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PREDICTION  OF  ATMOSPHERIC  MTF  AND  APPLICATION  TO  IMAGE 
RESTORATION,  BASED  ON  METEOROLOGICAL  DATA 

Fourth  Quarterly  Report 

The  Third  Quarterly  Report  contained  a  statistical  model  for  predicting  size 
distribution  of  coarse  aerosols  (radius  >0.16  iim)  according  to  weather  forecast.  This 
applies  to  arid  or  semi-arid  regions  in  which  the  soil  is  dry  so  that  soil-derived 
particulates  are  uplifted  by  the  wind.  The  smaller  fine  particulates  are  then  modelled 
with  MODTRAN  or  LOWTRAN  7. 

A  detailed  comparison  of  aerosol  MTF  prediction  with  aerosol  MTF 
measurement  indicates  the  reliability  of  this  model  is  excellent,  except  for  extreme  cases 
defined  in  The  Third  Quarterly  Report  in  which  the  aerosol  size  distribution  model  does 
not  hold.  Unfortunately,  quantitative  results  of  this  comparison  of  MTF  prediction 
with  measurement  are  not  yet  complete.  Qualitatively,  they  indicate  the  coarse  aerosol 
size  distribution  prediction  is  excellent. 

Here,  a  statistical  model  for  predicting  coarse  aerosol  size  distribution  according 
to  weather  for  areas  in  which  soil  moisture  prevents  much  uplifting  of  soil-derived 
particulates  is  presented  in  Appendix  1.  This  pertains  to  climates  which  are  not  arid  or 
semi-arid. 

Appendix  2  contains  a  step-by-step  procedure  to  calculate  aerosol  MTF,  based 

on  predicted  aerosol  size  distribution  (n(x)  in  step  10)  and  hardware  field-of-view  and 

dynamic  range. ^  Step  12  in  it  can  be  used  to  predict  atmospheric  transmission  ta 

according  to  weather  ,where 
L 

Xa=  exp[-  J  a[A-,n(x),z’]dz'.  (1) 

0 

In  (1),  a  is  medium  extinction  coefficient  calculated  in  step  12,  and  L  is  path  length. 
Thus,  the  aerosol  size  distribution  models  for  arid  and  non-arid  regions  can  be  used  in 
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Apependix  2  to  predict  both  aerosol  MTF  and  atmospheric  transmission  according  to 
meteorological  coditions. 

Finally,  the  image  restoration  work  has  continued,  and  a  comparison  of  various 
techniques  with  which  to  restore  images  degraded  by  atmospheric  blur  is  enclosed  as 
Appendix  3.  These  techniques  apply  to  visible  and  thermal  infrared  wavelengths,  and 
are  based  upon  predicted  or  measured  atmospheric  MTF.  The  improved  Wiener  filter 
presented  in  The  First  Quarterly  Report  yields  best  results  when  turbulence  blur  is 
significant.  When  aerosol  blur  is  more  dominant,  a  variety  of  techniques  all  yield 
similar  results.  In  all  cases,  the  inclusion  of  aerosol  blur  in  the  atmospheric  MTF  yield 
restoration  limited  by  hardware  alone,  as  if  there  were  no  atmospheric  blur. 

Reference 

1.  D.  Sadot  and  N.  S.  Kopeika,  "Imaging  through  the  atmosphere:  practical 
instrumentation-based  theory  and  verification  of  aerosol  MTF,"  J.  Opt.  Soc. 
Amer.  A,  vol.  10,  pp.  172-179,  Jan.  1993. 
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Statistical  model  for  aerosol  size  distribution  parameters  according  to  weather  parameters 
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ABSTRACT 

Predictions  of  atmospheric  transmittance  in  desert  aerosol  environments  using  MODTRAN  code 
diverge  significantly  from  measured  data.  Good  prediction  of  the  desert  particulate  size 
distribution  is  required  in  order  to  predict  atmospheric  scattering  and  absorption  parameters.  It  is 
also  essential  to  the  prediction  of  the  aerosol  atmospheric  modulation  transfer  function  which  is 
often  the  dominant  component  of  the  overall  atmospheric  MTF.*  Recently^  an  effort  to  predict 
statistics  but  not  size  distribution  according  to  simple  weather  parameters  has  been  made  for 
coarse  desert  aerosols.  A  quantitative  analysis  of  the  desert  particulate  size  distribution  models 
was  also  performed^.  In  this  research  the  size  distribution  parameters  measured  by  optical 
counters  are  related  to  weather  parameters.  Known  statistical  and  analytical  models  such  as 
MODTKAN  relate  the  size  distribution  parameters  only  to  relative  humidity  for  continental 
atmospheres.  Although  humidity  has  a  significant  role  in  the  prediction  of  aerosol  size  statistics, 
other  weather  parameters  are  seen  here  to  strongly  influence  also  the  size  distribution  parameters. 

Comparisons  such  as  the  above  can  be  used  to  predict  under  which  conditions  the  MODTRAN 
aerosol  models  have  good  or  poor  accuracy.  It  is  also  hoped  that  they  will  lead  to  improvements 
in  MODTRAN,  improving  the  accuracy  of  the  humidity  dependence  as  well  as  by  incorporating 
other  meteorological  parameters  into  the  MODTRAN  prediction  models. 


1.  Introduction 


The  aerosols  that  are  suspended  in  the  atmosphere  are  a  major  factor  that  influence  the  quality  of 
images  that  are  propagating  through  the  atmosphere.  Scattering  and  absorption  of  light  rays  by 
the  suspended  aerosols  contribute  to  the  loss  of  contrast  of  the  received  image,  forward  light 
scattering  in  the  direction  of  the  image  propagation  also  reduce  the  atmospheric  modulation 
transfer  function  (MTF)  and  cause  to  image  blur. 

Prediction  of  the  aerosol  size  distribution  can  lead  to  the  prediction  of  the  scattering  and 
absorption  coefficients  and  also  the  scattering  phase  function  of  the  atmospheric  propagating 
medium.  These  properties  characterize  the  extinction  coefficient  of  the  atmosphere  as  well  as  the 
MTF  of  the  atmosphere. 

2.  Experiment 

The  particulate  size  distribution  of  the  atmosphere  was  measured  over  the  radius  range  of  0.16 
to  10  |im  using  a  PMS  CSASP-100  optical  particles  counter.  The  optical  counter  was  located  at  a 
height  of  25  m  above  ground.  Integration  time  of  half  an  hour  was  chosen  and  every  half  an  hour 
the  average  particulate  size  distribution  was  recorded.  Meteorological  parameters  were  ^sp 
measured  close  to  the  PMS  location  and  the  half  hour  average  of  the  temperature,  relative 
humidity,  wind  speed,  wind  direction  and  solar  flux  were  also  recorded.  The  measurements  were 
taken  at  the  Ben-Gurion  University  of  the  Negev  in  Beer-Sheva,  Israel.  Data  was  collected  form 
June  1991.  Both  PMS  counter  and  Campbell  Scientific  meteorological  station  were  connected  to 
a  personal  computer. 

Since  this  research  is  trying  to  relate  the  aerosol  particulate  size  distribution  parameters  to 
simple  weather  parameters,  the  integration  time  selected  for  counting  and  classifying  particles 
according  to  their  size  was  to  be  shorter  than  the  meteorological  time  constants.  On  the  other 
hand,  a  short  integration  time  leads  to  a  particulate  size  distribution  curve  with  insufficient  data 
which  is  hard  to  fit  to  a  known  size  distribution  model.  In  our  measurements  an  integration  time 
of  half  an  hour  was  chosen  as  a  compromise  between  the  need  for  a  smooth  size  distribution  curve 
and  the  need  to  measure  the  aerosol  statistics  under  stationary  meteorological  conditions. 

Fig.  1.  describes  the  changes  in  the  standard  deviation  of  the  temperature  during  a  typical  summer 
day.  The  standard  deviation  here  is  calculated  for  measurement  time  periods  of  half  an  hour. 
Fig.  2  shows  the  changes  in  the  standard  deviation  of  the  relative  humidity  during  the  same  day. 
These  figures  show  that  the  most  rapid  changes  occur  after  sunrise  and  sunset.  In  our 
measurements  the  standard  deviation  of  the  temperature  did  not  exceeded  2°C  and  that  of  the 
relative  humidity  did  not  exceed  5  %.  In  most  of  the  cases  the  standard  deviations  were  much 
more  less  than  these  values. 
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Fig.  1.  Half  hour  standard  deviation  of  air  temperature  as  a  function  of  time  during  a  typical 
summer  day. 


Fig.  2.  Half  hour  standard  deviation  of  relative  humidity  as  a  function  of  time  during  a  typical 
summer  day. 

3.  Data  Analysis 

The  deviations  in  the  aerosol  size  distribution  are  large.  The  aerosols  in  the  area  of  our 
measurements  consist  of  a  large  variety  of  aerosol  sources.  Aerosol  chemistry  and  area 
meteorology  are  described  in  detail  in  ref.  2.  The  aerosol  sources  are  urban  and  rural,  marine, 
haze,  industrial  and  desert  dust  aerosols  that  are  transferred  from  the  Sahara  desert  and  from  the 
Arabian  desert^  during  dust  storms, 
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Fig.  3.  describes  two  extreme  aerosol  size  distributions  that  were  measured  during  the  year 
1992.  The  differences  between  the  two  extreme  cases  is  about  two  orders  of  magnitude.  The 
lower  curve  is  not  continuous  because  particles  at  several  radii  were  not  counted  during  an 
integration  period  of  half  an  hour.  Fig.  4.  shows  the  normalized  standard  deviation  of  the 
measured  aerosol  at  each  radius.  The  standard  deviation  ranges  from  0.5  at  radius  of  0.16  |im  to 
12  at  radius  of  10  |im.  The  aerosol  size  distribution  consists  of  several  distributions  that  can  be 
modeled  by  superposition  of  lognormal  distributions.^  The  standard  deviation  of  the  particle  count 
is  smaller  for  small  radii  than  for  high  radii  since  the  low,  or  fine  mode  of  the  aerosol  distribution 
consists  mainly  of  local  aerosols  while  the  large  aerosol  population  contains  also  particles  that  are 
transported  from  remote  areas.  The  local  aerosols  are  less  affected  by  changes  in  the  weather 
parameters  while  the  larger  particles  are  very  sensitive  to  the  wind  speed  and  to  the  wind  direction 
which,  in  some  cases,  yields  information  concerning  the  source  of  the  remote  aerosols. 


A  simple  size  distribution  that  can  characterize  the  coarse  aerosol  size  distribution  is  the  Junge 
power  law  distribution.^  The  size  distribution  is  defined  by  only  two  parameters  A  and  a: 

A  ~0-  / 1  \ 

—  =  Ar  (1) 

dr 

r  is  the  particle  radius  and  n  is  the  particulate  number  concentration. 

The  advantage  of  this  size  distribution  is  that  a  measured  size  distribution  can  be  fitted  to  this 
model  by  a  simple  one  variable  linear  regression. 


In  our  analysis  the  measured  size  distributions  were  fitted  to  Junge  power  law  models.  The  data 
were  separated  into  two  major  groups.  One  consists  of  the  measured  data  during  the  summer 
season  from  May  to  October.  This  corresponds  to  arid  and  semi-arid  regions.  The  second  group 
consisted  of  data  that  was  measured  during  the  winter  season,  starting  from  November  until  April. 
This  corresponds  to  wet,  non-arid  regions. 

The  climate  in  our  region  is  such  that  there  are  two  major  seasons.  Summer  is  the  hot  and  dry 
season  with  no  participation  at  all.  Winter  is  the  colder  season.  All  participation  is  during  this 
season. 


3.1  A  model  for  the  summer  season 

Fig.  5.  describe  the  average  particulate  size  distribution  and  the  power  law  fit  All  6000 
measurements  of  the  aerosol  size  distribution  from  May  1992  to  October  1992  were  fitted  to  the 
power  law  distribution  and  the  values  of  a  and  A  were  extracted.  The  average  square  of  the 
correlation  coefficient  was  0.95. 

The  power  law  size  distribution  parameters  A  and  a  were  related  to  the  weather  parameters  that 
were  measured  nearby  the  particulate  counter  using  multiple  linear  regression  analysis.  The 
summer  model  obtained  for  A  is  given  by: 

A  =  aQ-k-a^rh+a2rl^ +a.^r}^ +a^rh*  (2) 

were  rh  is  the  relative  humidity.  The  coefficient  of  the  model  are  given  in  table  1. 
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Table  1.  Summer  model  for  the  prediction  of  the  parameter  A. 


an 

ai 

_ 

a. 

a4 

1.917 

-0.128 

4.61- 10-3 

6.85- 10-5 

3.73-10-7 

The  model  obtained  for  a  was  the  following; 

(X=bQ+b^rh+b2rh+b3WS  (3) 

were  ws  is  the  average  wind  speed  in  m/s.  The  coefficients  for  the  prediction  of  a  are  given  in 
table  2. 


Table  2.  Summer  model  for  the  prediction  of  the  parameter  a. 


bn 

b, 

^2 

b. 

2.68 

2.3-10-3 

8.16-10-5 

-2.65-10-2 

The  A  parameter  in  this  model  depends  only  on  the  relative  humidity.  Other  dependencies  were 
weak  and  insignificant.  The  parameter  a  depends  also  on  the  average  wind  speed.  Fig.  6  describes 
in  an  ascending  order  the  measured  a  and  the  modeled  a.  From  this  figure  it  is  shown  that  there  is 
a  good  agreement  between  the  measured  parameter  and  the  modeled  parameter  except  for 
extreme  weather  conditions.  Fig.  7.  shows  the  measured  parameter  A  in  ascending  order  and  the 
modeled  parameter.  The  model  gives  good  correlation  here  too  except  for  extreme  weather 
conditions. 

Figs.  8.  show  typical  measured  particulate  aerosol  distributions  (solid  lines)  and  predicted  size 
distribution  according  to  the  measured  weather  parameters  (dashed  lines).  The  curves  shows  that 
there  is  a  very  good  agreement  between  the  measured  data  and  the  model  for  the  parameters  a 
and  A  under  a  variety  of  weather  conditions. 


Fig.  3.  Average  aerosol  size  distribution  curve  and  size  distributions  with  the  property  of 
minimum  and  maximum  number  of  particles  that  were  measured. 
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Fig.  4.  Normalized  standard  deviation  of  the  aerosol  size  distribution. 


Fig.  5.  The  measured  average  particulate  size  distribution  and  its  power  law  fit. 
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Fig.  6.  Measured  a  parameter  and  predicted  a. 


Fig.  7.  Measured  A  parameter  and  predicted  A. 
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Fig.  8a.  Measured  particulate  size  distribution  and  the  modeled  size  distribution. 


Fig.  8b.  Measured  particulate  size  distribution  and  the  modeled  size  distribution. 


Fig.  8c.  Measured  particulate  size  distribution  and  the  modeled  size  distribution. 

4.2  A  model  for  the  winter  season 

A  similar  approach  for  modeling  the  aerosol  size  distribution  for  the  winter  season  which  is  the 
wet  and  rainy  season  in  our  area  was  used.  The  two  parameters  of  the  power  law  model  for  the 
aerosol  size  distribution  were  related  to  weather  parameters.  The  major  difference  between  the 
two  seasonal  models  is  that  for  the  winter  season,  not  only  were  the  measured  weather  parameters 
necessary  but  also  the  history  of  some  of  them.  The  model  here  is  based  on  about  5000  different 
measurements  of  the  aerosol  size  distribution  during  winter  seasons. 

The  model  for  In  A  is  the  following: 

In  A  =  flg  4-aj  rh+Uj  rh2i  +a^  {rhiA  +a^  w.s+  a^wsiiA  (iwm  +a^  tm+a^  sf  (4) 

where  rh  is  relative  humidity,  r/124  is  the  last  24  hour  average  of  the  relative  humidity,  ws  is  the 
wind  speed  in  m-s'^  m24  is  the  24  average  of  the  wind  speed,  tm  is  the  temperature  in  ®  C,  and  sf 
is  the  solar  flux  in  kwatt-m'^. 

The  following  model  was  obtained  for  a: 

a  =  Z?o +!>j  r/i+^2  ^^24  +b^{rh24  +b^ws+b^ws2A  +b^(ws2A)^ +bj  tm+b^sf  (5) 

Fig.  9  show  an  example  of  the  measured  aerosol  size  distribution  and  the  modeled  size  distribution 
according  to  the  weather  parameters  and  the  24  hour  averages  of  the  wind  speed  and  the  relative 
humidity. 
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Fig.  9.  Measured  particulate  size  distribution  and  the  modeled  size  distribution  for  the  winter 
season. 
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5.  Conclusions 


The  major  conclusion  of  this  research  is  that  the  coarse  aerosol  size  distribution  parameters  can  be 
predicted  according  to  simple  meteorological  parameters.  A  strong  statistical  relation  exists 
between  the  aerosol  size  distribution  curve  and  the  local  weather.  Although  the  Junge  power  law 
distribution  is  not  accurate  since  it  does  not  model  the  size  distribution  at  radii  that  are  smaller 
than  0.1  p,m  and  the  fine  details  of  the  curve  are  lost,  it  can  model  roughly  the  aerosol  size 
distribution  and  it  is  very  easy  to  extract  its  two  parameters  A  and  a  even  if  the  aerosol  histogram 
is  composed  of  relatively  low  counts  of  particles.  The  strong  statistical  dependence  of  the  power 
law  distribution  on  local  weather  parameters  show  us  that  in  most  cases  the  statistics  of  the 
aerosols  are  influenced  by  the  local  weather.  There  are  few  exceptions,  mainly  when  the  aerosol 
size  distribution  consists  of  particles  such  as  desert  dust  particles  or  marine  particles  that  are 
transported  from  other  areas. 

Introducing  the  history  of  the  meteorological  parameters  such  as  the  previous  24  hours  average  of 
the  wind  speed  and  that  of  the  relative  humidity  was  very  significant  to  the  winter  model  and 
contributed  significantly  to  the  correlation  of  the  model  with  measurement. 
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Appendix  2 


Calculating  Practical  Aerosol  MTF 

Mie  scattering,  absorption  and  extinction  coefficients 
Mie  phase  function 


1.  Define  the  following  parameters 

m  -  relative  refractive  index  (complex). 
X  -  radiation  wavelength  (meters). 
a  -  particle  radius  (meters). 

Ik  a 

X  =  — - size  parameter. 

p  =  x-rn. 

6  -  scattering  angle  (steradians). 

|i  =  cos6. 

n{x)  =  particle  size  distribution 


2.  Create  the  spherical  Bessel  functions,  which  satisfy  the  recurrence  relations 

Zn-l(p)  +  Z/!+l(p)  = - t:n(p).  (1) 

P 

^w(p)  ~  Zn— 1  (p)  ~  l)  Zn+1  (p)  ’  (2) 

dp 

where  is  either or  y„  From  the  first  two  orders 


./X  sinp  COSO 

p^-p- 

(3) 

>o(p)=-“^P. 

^^(p)_cosp_sinp_ 

(4) 

hiffher-order  functions  can  he  venerated  hv  recurrence. 

3.  Create  the  spherical  Bessel  functions  of  the  third  kind  (Hankel  functions) 

hnip)  ~  j  nip)  ~  ^  y 


(5) 


4.  Create  the  Riccati-Bessel  functions 

<p„(p)  =  p;„(p),  ^„(p)  =  p/zn(p).  (6) 

5.  Calculate  the  scattering  parameters  a„  and 

_»?(pn(p)(plU)-(pn(-y)(pl(p) 

w9„(p)CW-^«(^)9„(p)  ’  ^ 

_  9n(p)9„(j^)-»^(P„(-v)<pi,(p) 


6.  Calculate  the  Mie  scattering  ,  extinction  and  absorption  efficiency  factors 

Q^c{^n,x)  =  —  E(2«-f  l)(|a„p  +  |^,ip),  (9) 

X  n=l 

Qext  =  —  x(2n  + 1)  Re{a,.  +  ^„},  (10) 

A'  n=l 

Q„bs=Qe.t-Qsc-  (H) 


7.  Create  the  angle-dependent  functions  which  satisfy  the  recurrence  relations 

_  _  2/7  -  1  _  n 

ttn  ~  ~Kn-2^ 

n-l  /7  -1 

tn  ~  ^l-l  Kn  ~  l) 

where  71q  =  0  and  k  i  =  1. 


(12) 

(13) 


8.  Calculate  the  complex  amplitudes  of  scattered  electric  fields 


5i(/77,A,e)=  i  .AantZn+bnXn), 

n=in\n  +  Ij 

(14) 

6’2(«r,-^,6)  =  X  /  ^\{bnTtn  + anX^ 

„=in(/7  +  lj 

(15) 

where  Si  and  Sz  are  components  perpendicular  and  parallel  to  the  scatter  plane. 


9.  Calculate  the  dimensionless  intensity  parameters  (for  scattered  light  only) 

il(A,/77,e)  =  5i5p 

i2[x,m,^)  =  S2S2- 


(16) 

(17) 


10.  Use  the  measured  atmospheric  particle  size  distribution  to  calculate  the  medium's  scattering 
coefficient  (scattering  cross  section  per  unit  volume) 


(18) 


=  T^k  ^Iq  x^n{x)Q^^{x)dx. 

11.  Calculate  the  medium's  phase  function 

^(®)=-J^F”(^)[n(®)+/2(0)]^-«.  (19) 

which  obeys  the  normalization  condition: 

P(e)rAv  =  47t  (20) 

where  the  integration  is  with  respect  to  solid  angle  over  all  directions  corresponding  to  all 
space  around  a  point. 

12.  Calculate  the  medium's  extinction  coefficient  (extinction  cross  section  per  unit  volume) 

<3[k,n{x)]=^TZk~^\Q  x'^n{x)Q^^,{x)dx,  (21) 

13.  Calculate  the  medium's  absorption  coefficient  (absorption  cross  section  per  unit  volume  ) 

Aa[X,/7(x)]  =  x'^n{x)Qabs{x)dx,  (22) 

Atmospheric  absorption  -  spatial  frequency  (dependence 


fa  <  feu 
fa  >  fac. 

(23) 

15.  Calculate  the  new  MTF  asymptote  at  high  spatial  frequencies  =  Ma  (~)  (differing  from  the 
atmospheric  transmittance),  and  calculate  the  effective  optical  depth 

r.j=-ln[Md4]-  (24) 

Practical  aerosol  MTF  ' 

16.  Approximate  the  medium's  phase  function  by  a  Gaussian  form  (not  necessary)  ■ 

P(e)  =  4.^/J^exp{-ape2}.  (25) 

17.  Calculate  the  "classical"  aerosol  MTF  using  the  effective  optical  depth 


14.  Calculate  aerosol  MTF  including  absorption  spatial  frequency  dependence 


^«(/a)  =  exp<  -  5az[-^  I  >-exp^  exp^-i'^c  1- 
W  ac)  V 

exp{-  5«s}-exp{[l-  exp{-  5nz}](-  A«z)}, 


(26) 


(f(/ J  = 


exp 


0 


1  -  exp 


J  a 


a 


P  J 


\dz' 


18.  Determine  the  main  instnimentation  limitations: 

a  -  Finite  angular  spatial  frequency  bandwidth  -  /^niax- 
b  -  Finite  field  of  view  -  ©'max- 

c  -  Dynamic  range  limitation  -  minimal  detected  received  irradiance,  which  is  related  to  maximal 
angle  of  detected  scattered  light  -  ©"max- 

d  -  SNR  limitation  -  minimal  detected  received  irradiance,  which  is  of  higher  intensity'  than  the 
background  noise  -  ©'"max- 


19.  Choose  the  worst  case  between  (18b)  -  (18d),  obtaining  the  instrumentation  limitation  of  maximal 
angle  of  recorded  scattered  light  -  ©^ax  =  niin  {©'niax-©"max.©"'max} 

20.  Calculate  the  "classical"  specific  intensity  function  for  the  case  of  an  infinite  plane-wave 

/  \  1  ^  /  \  /  \ 

— o  i  df^-exp(-jefJ-K(fJ.  (27) 


20a.  Iri  case  where  the  system’s  MTF  is  known,  (27)  can  be  replaced  by 

‘'(S)  =  7-^  (27a) 

\2n)  -oc 

where  Moifa)  system  MTF. 


2 1 .  Calculate  the  practical  aerosol  MTF  by  inverse  Fourier  transforming  the  specific  intensity  function 

PMTF„{fJ=  ®r'(e'«p(;e/J/(e). 

“  ©max 


(28) 
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Comparison  between  high-resolution  restoration  techniques 
of  atmospherically  distorted  images 
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Department  of  Electrical  and  Computer  Engineering 
Ben-Gurion  University  of  the  Negev 
Beer-Sheva,  Israel 

Abstract 

A  system  approach  is  applied  to  overcome  atmospheric  degradation  of  remotely  sensed 
images.  A  comparison  is  presented  between  different  filtering  techniques  for  restoration  of 
distorted  images  for  both  the  visible  and  thermal  infrared  spectral  regions.  Restoration  methods 
include  spatial  and  spatial  frequency  filters.  Best  results  are  obtained  by  using  a  fractal  model 
to  describe  the  image's  power  spectral  density  for  scenes  in  the  visible  spectral  range. 
Atmospheric  effects  are  best  modeled  by  a  noisy  spatial  frequency  filter  composed  of  an 
average  component  described  by  the  average  atmospheric  and  hardware  modulation  transfer 
function,  and  a  noisy  component  modeled  by  the  atmospheric  point  spread  function's  power 
spectral  density.  The  most  impressive  restorations  in  the  visible  range  are  achieved  by 
combining  between  the  last  and  the  fractal  model  for  the  object's  power  spectral  density.  In  the 
thermal  range,  however,  several  restoration  techniques  yielded  very  good  results,  with  no  one 
technique  most  advantageous  over  the  others.  The  methods  presented  here  are  capable  of 
yielding  real-time  image  restoration  with  resolution  limited  essentially  only  by  the  hardware 
itself  in  both  wavelength  regions. 
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1.  Introduction 

Imaging  through  the  atmosphere  has  progressed  significantly  in  the  last  decade,  in  both 
the  visible  and  infrared  spectral  ranges.  It  is  the  atmosphere  which  usually  limits  image  quality, 
particularly  for  long  atmospheric  paths.  The  main  atmospheric  distortions  are  caused  by  optical 
turbulence,  and  scattering  and  absorption  by  particulates  in  the  atmosphere.  The  main  effect  of 
the  turbulent  medium  over  long  exposures  is  to  produce  wavefront  tilt,  which  causes  image 
shifts  at  the  image  plane.  These  can  be  whole  image  shifts,  or  different  shifts  of  different  parts 
of  the  image,  depending  on  the  isoplanatic  patch,  according  to  the  turbulence  parameters 
(turbulence  strength,  inner  and  outer  scales).  The  image  distortions  caused  by  the  wavefront 
tilt  (typically  on  the  order  of  tens  or  maybe  hundreds  of  microradians)  can  be  partially 
compensated  for  either  by  adaptive  optics  techniques  or  by  using  sufficiently  short  exposure 
time,  less  than  the  characteristic  fluctuation  time  (usually  a  few  milliseconds).  Another  way  is 
by  using  a  wavefront  sensor  to  construct  the  image  phase  in  addition  to  image  intensity,  and 
therefore  restore  the  image  in  a  deterministic  way^.  This  is,  however,  not  practical  for  real  time 
restoration.  Turbulence  distortion  effects  are  characterized  in  both  short  and  long  exposures  by 
modulation  transfer  functions  (MTFs)2*3.  These  are  statistical  averages  of  turbulence  random 
processes  which  have,  in  addition,  a  nonnegUgible  short-time  variance  which  affects 
significantly  image  blur.  In  addition  to  turbulence,  there  are  scattering  and  absorption  effects 
produced  by  molecules  and  aerosols  in  the  atmosphere.  These  cause  both  attenuation  and 
image  blur,  according  to  the  atmospheric  aerosol  MTF^-^.  Unlike  turbulence,  the  aerosol  MTF 
affecting  the  image  is  also  related  strongly  to  instrumentation  limitations^-®-^  on  the  angles  of 
light  scatter  actually  recorded  in  the  image.  Thus,  aerosol  effects  can  be  approximated  by  an 
MTF  of  limited  long-time  variation,  according  to  weather  conditions. 

Here,  image  processing  techniques  are  aided  by  the  knowledge  of  the  corruption  from 
which  the  image  has  suffered  while  propagating  through  the  atmosphere.  These  techniques  are 
preferable  rather  than  standard  blind-deconvolution  techniques,  as  long  as  the  corrupting 
medium,  i.e.  atmosphere,  is  modeled  sufficiently  accurately.  One  way  of  modeling  the 


3 


atmospheric  distortions  is  by  the  average  atmospheric  MTF.  In  this  technique,  atmospheric 
effects  are  modeled  by  a  deterministic  filter.  This  method  yields  good  restoration  results  when 
atmospheric  MTF  does  not  vary  too  much,  particularly  when  turbulence  effects  are  not 
dominant.  This  is  the  case  in  the  thermal  range,  or  while  only  low  spatial  frequencies  are 
involved.  When  turbulence  is  dominant,  however,  modeling  atmospheric  MTF  by  a 
deterministic  filter  is  insufficient,  due  to  significant  MTF  variance.  Atmospheric  effects  in  that 
case  should  be  modeled  by  a  noisy  spatial  frequency  filter  including  an  average  component 
described  by  the  average  atmospheric  MTF,  and  a  noisy  component  which  contains  the  MTF 
variance. 

The  power  spectral  density  (PSD)  of  the  original  picture  is  often  required  in  order  to 
apply  the  restoration  filter  to  the  corrupted  image  (as  in  the  standard  Wiener  filter,  for 
example).  One  way  of  estimating  this  term  is  by  rigorously  assuming  that  it  equals  the  image's 
PSD.  This  is  of  course  not  true,  but  in  many  cases  this  is  a  sufficiently  good  assumption.  To 
check  this  issue,  one  can  restore  the  image  iteratively.  This  process  should  stop  when  restored 
images  converge.  Usually,  it  happens  after  one  or  two  iterations.  Another  way  of  estimating 
the  object's  PSD  is  by  assuming  it  obeys  a  fractal  model.  It  wOl  be  shown  here  that  using  a 
fractal  estimation  model  yields  very  good  results  for  visual  images. 

In  this  paper,  the  fractal  dimension  estimation  process  is  presented.  Then,  the  improved 
Wiener  filter  technique  is  introduced.  This  is  followed  with  other  spatial  filtering  restoration 
techniques.  Results  of  the  different  methods  wiU  be  compared,  including  some  examples  of 
restored  thermal  images  as  well  as  visible  images. 

2.  Fractal  model 

If  it  were  possible  to  know  characteristics  of  the  original  background  with  which  one 
was  dealing,  one  could  potentially  design  more  efficient  and  effective  algorithms  to  perform  the 
restoration.  In  this  paper,  we  describe  the  restoration  of  images  for  natural  scenes.  It  has  been 
proposed  that  natural  scenes  are  fractal  in  nature;  pictures  taken  in  the  visual  spectral  range 
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will  also  be  fractal^o  contrast  to  infrared  images  which  tend  to  be  modeled  by  Markov 
random  fields 

Our  restoration  of  the  scenes  based  on  their  similarity  to  fractal  models  is  performed  as 
follows  Fractal  scenes  in  the  visual  regime  have  power  spectra  with  a  —  frequency 

dependence.  We  assume  that  even  in  real  natural  pictures,  the  behavior  will  be  correct  for 

wide  bands  of  frequencies.  The  types  of  blurring  which  occur  will  affect  the  higher  frequencies 
much  more  than  the  lower  frequencies.  We  rather  arbitrarily  assume  that  the  corruption  has 
only  occurred  to  the  highest  frequencies  for  which  0.1%  of  the  energy  of  the  picture  is  found; 

although  this  seems  minuscule,  in  a  large  number  of  picture  distortions  this  is  indeed  where  the 

major  chanses  occur.  The  remaining  power  of  the  lower  frequencies  are  fit  to  —  where  n  is 

/" 


determined  by  a  least-square  fit  of  the  log  of  the  radial  power  to  the  frequency.  We  then 
assume  that  the  higher  frequencies  in  the  original  picture  originally  followed  the  same  — 

dependence  before  corruption.  These  frequencies  are  then  restored. 

Results  of  images  restored  by  using  best  fits  between  the  object's  PSD  to  fractal  model 
are  presented  in  Fig.  1.  It  is  clear  that  scene's  resolution  is  improved.  This  method  yields  very 
good  results  for  pictures  which  include  relatively  "natural"  scenes.  However,  when  the  scenes 
become  more  "man-made",  restoration  efficiency  decreases.  In  the  thermal  range,  however,  the 
fractal  model  completely  failed,  as  demonstrated  in  Fig.  2.  Both  "man-made"  and  thermal 
images^  ^  are  known  not  to  be  fractal  in  nature,  but  rather  Markovian;  thus  failure  of  a  fractal 
method  is  not  unexpected. 


3.  Improved  Wiener  filter 


The  improved  Wiener  fUter^^  is  based  on  the  standard  Wiener  filter^^^  ^n 
additional  noise  model  deriving  from  the  atmospheric  MTF  variance.  As  explained  in  the 
introduction,  the  atmospheric  MTF  derives  from  both  turbulence  and  aerosol  forward 
scattering  effects.  While  aerosol  MTF  is  fairly  constant  as  long  as  atmospheric  conditions  do 
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not  vary  too  much,  turbulence  MTF  changes  with  time  due  to  it's  tilt  jitter  characteristic.  These 
tilts  are  random  and  their  temporal  power  spectra  are  usually  limited  to  several  tens  up  to  a  few 
hundred  hertz  under  ordinary  atmospheric  conditions.  The  image  distortions  caused  by 
atmospheric  MTF  are  thus  regarded  as  the  sum  of  a  deterministic  and  a  random  filteri^-i^.  The 
deterministic  filter  includes  aerosol'^-^  and  average  turbulence  MTFs^-^,  while  the  random  filter 
includes  the  noise  component  induced  on  the  imaging  system,  both  by  the  turbulence  MTF 
variance  and  hardware^^.  Stated  mathematically,  the  atmospheric  filter  is  determined  by 

h'=h  +  ni  (1) 

where  /z'  is  the  instantaneous  atmospheric  point  spread  function  (PSF),  h  is  the  average 
atmospheric  PSF,  and  n/  is  an  additive  random  component  with  zero  expectation.  Using  this 
model,  the  image  received  at  the  imaging  system  after  propagating  through  the  atmosphere  is: 

g(x,y)={h(x,y)  +  ni(x,y))®  f(x,y)  +  n2(x,y)  (2) 

where  g  is  the  received  image, /is  the  object,  x  and  y  are  transverse  spatial  coordinates,  and  «2 
is  an  additive  noise  imposed  by  the  instrumentation,  including  optics,  digitization,  electronics 
etc.,  but  not  by  the  atmosphere.  Fourier  transforming  (2)  yields: 

G(u,v)=[H(u,v)  +  Ni(u,v)]-  F(u,v)  +  N2(u,v)  (3) 

where  G,  H,  F,  Nj  and  N2  are  Fourier  transforms  of  g,  h,f,nj  and  ^2  respectively,  and  u  and  v 
are  spatial  frequency  coordinates.  The  received  image  thus  is  a  sum  of  a  deterministic  part  Gj 
and  a  random  part  N 

G  =  G\  +  N  (4) 

where  G\{u,v)  =  h{u,v)- F{u,v)  (5) 

and  /^(m,v)=  F(m,v)- Ari(M,v)  + A/'2(m.v).  (6) 

The  improved  Wiener  filter  is  defined  similarly  to  the  standard  Wiener  filter,  differing  by  the 
noise  component  which  includes  an  additional  term  imposed  by  the  atmosphere^^ 

_ \H{u,vf\ _ 

H («,v)  •(li?(M, v)^I+[s„,„,(m,v)  +  5„2„3(m,v)  /  S ff  (m,v)]) 


(7) 
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where  M  is  the  restoring  filter,  H{u,v)  is  the  average  atmospheric  MTF, 

Sff{u,v)  are  the  power  spectral  densities  of/,  n\,  and  ni,  where  n\  and  ni 
are  the  inverse  Fourier  transforms  of  Ni  and  Ni  iii  (6). 

Assuming  independence  between  aerosol  and  turbulence  effects,  the  term  h{u,v)  can  be 

measured  or  calculated  by  a  multiplication  of  the  turbulence  MTF  (either  short  or  long 

exposure  case)  and  aerosol  MTF.  Turbulence  MTF  can  be  evaluated  with  the  knowledge  of 
standard  meteorological  parameters  using  a  prediction  modeF^  (verified  independently  by 

U.S.  Army  Night  Vision  Laboratory)  or  IMTURB  or  PROTURB  (U.S.  Army  Atmospheric 
Sciences  Lab.),  and  the  aerosol  MTF^  can  be  evaluated  according  to  knowledge  of  particle  size 

distribution,  which  can  also  be  predicted  via  LOWTRAN,  MODTRAN,  or  other  models^*.  The 
term  Sff{u,v)  can  be  estimated  either  by  using  the  received  image  G{u,v)  or  by  using 

estimation  models  of  the  object’s  PSD,  which  has  been  shown  to  obey  a  fractal  model  in  the 
visible  range,  and  a  Markovian  model  in  the  thermal  range^k  The  term  is  assumed 

to  be  constant  for  aU  spatial  frequencies  since  the  additive  noise  n2  is  assumed  to  be  white 

noise.  This  assumption  is  commonly  used  and  very  practical,  and  it  has  a  relatively  weak  effect 
on  the  Wiener  filter.  The  term  v)  is  very  important,  since  it  includes  the  random  part  of 

the  atmospheric  distortions.  One  way  of  estimating  is  by  a  direct  measurement.  By 

using  the  relation 

S^^{u,v)=^E[n}{u,v)},  (8) 

and  the  Fourier  transform  of  (1),  i.e., 

Ni{u,v)  =  H'{u,v)-H{u,v),  (9) 

it  follows  that  equals  the  variance  of  the  instantaneous  atmospheric  MTF: 

5,^^(m,v)  =  £{H’2(m,v)}-h2(m,v).  .  (10) 

Since  the  contribution  to  the  random  part  of  the  atmospheric  MTF  is  due  mainly  to  turbulence 
rather  than  aerosols,  in  (8)-(10)  atmospheric  MTF  refers  to  turbulence  only.  However,  in  (7), 
H{u,v)  includes  aerosol  MTF  in  addition  to  turbulence  MTF  since  it  refers  to  the  average 
atmospheric  MTF.  The  variance  of  H'  can  be  evaluated  by  calculating  both  terms  of  the  right 
hand  side  of  (8).  This  can  be  carried  out  by  measuring  a  series  of  instantaneous  atmospheric 
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MTFs,  and  evaluating  the  average  of  both  the  MTF  and  its  square.  This  is,  however,  not  a  very 
practical  way,  particularly  when  real  time  image  restoration  is  concerned. 

An  alternative  way  of  obtaining  5„,„,(m,v)  is  by  evaluating  both  terms  of  the  right  hand 

side  of  (10).  The  second  term  IP  is  the  square  of  the  turbulence  MTF  which  can  be  predicted!^ 
or  measured.  The  first  term  e[h'^]  was  evaluated  analytically^^  and  yields: 

e{h'^{u,v)]  «  (M,v)Jj  {u' ,v')H{u'+u,v'+v)H(u'-u,v'-v)du'dv\  (11) 

oo 

Equation  (11)  determines  the  expected  value  of  the  squared  MTF,  or  in  other  words  the  Point 
Spread  Function's  Power  Spectral  Density.  The  integral  in  (11)  can  be  evaluated  numerically 
with  the  use  of  the  average  turbulence  MTF  only. 

Examples  of  restored  images  using  this  technique,  obtained  over  a  6.5  km  horizontal 

path  length  through  the  atmosphere,  are  presented  in  Fig.  3.  The  restoration  was  carried  out  by 
the  use  of  (7).  The  term  h{u,v)  was  the  measured  atmospheric  MTF.  The  term  S ff  {u,v)  was 

estimated  by  best  fit  to  a  fractal  model.  Obtaining  the  term  S ff{u,v)  from  the  received  image 
iteratevely  yielded  less  impressive  results.  The  term  assumed  to  be  white  noise 

and  the  term  was  evaluated  via  (11).  Restoration  time  was  only  about  2  seconds  per 

frame.  This  can  be  shortened  to  a  fraction  of  a  second,  using  parallel  processing  techniques 
already  available.  Therefore  restoration  via  this  method  can  be  in  real  time.  There  is  a  distinct 
improvement  in  fine  details  of  the  images,  even  though  the  image's  SNR  is  not  degraded 
significantly.  This  is  so  in  spite  of  the  severe  imaging  conditions  (long  horizontal  distance), 
where  turbulence  isoplanatic  patch  was  less  than  image  size,  and  standard  Wiener  filters  failed 
completely  in  trying  to  restore  the  image,  as  can  be  seen  in  Fig.  4.  In  the  standard  Wiener  filter, 
refers  to  white  noise  only. 

The  improved  Wiener  technique  yielded  very  good  results  in  the  thermal  range  too  (8- 
12  |im  spectral  window),  and  restoration^  for  a  picture  obtained  over  a  2  km  horizontal  path 

length  through  the  atmosphere,  using  this  technique,  are  presented  in  Fig.  5.  The  term 
Sff{u,v)  was  estimated  from  the  received  image,  and  no  iteration  was  needed  (restoration 

converged  after  the  first  iteration).  Since  the  thermal  images  involved  significantly  lower 
resolution  than  those  in  the  visible  range  (due  to  instrumentation  limitations),  the  main 
atmospheric  distortions  were  due  to  aerosol  rather  than  turbulence  MTF,  Therefore,  as 
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explained  previously,  MTF  variance  was  significantly  weaker  and  restoration  was  possible 
9 

using  an  average  atmospheric  MTF  (rather  than  instantaneous),  as  demonstrated  at  the  bottom 
right  of  Fig.  5.  This  is  important  while  considering  real  time  restoration  since  atmospheric 
MTF  can  be  estimated  once  to  create  an  efficient  restoration  filter  which  can  be  used  for  a 
large  set  of  images.  The  fact  that  atmospheric  MTF  in  the  thermal  range  is  relatively  static 
causes  relatively  impressive  restoration  results  using  standard  Wiener  filters  (which  failed  in  the 
visible  range).  An  example  can  be  seen  in  Fig.  6. 

4.  Constrained  filter 


The  constrained  filter^^  uses  a  prior  knowledge  about  the  picture  corruption  i.e. 
atmospheric  distortions,  in  addition  to  many  kinds  of  constraining  conditions,  in  order  to 
achieve  best  restoration  performance.  The  constraining  conditions  can  also  depend  upon  prior 
knowledge  about  the  target.  In  that  case,  improved  restoration  results  are  expected,  but  the 
filter  must  be  adapted  to  the  kind  of  scene  which  is  to  be  restored.  Here,  a  more  general 
condition  is  used,  which  is  restoration  using  the  atmospheric  MTF,  under  the  condition  of 
maximal  smoothness.  This  condition  minimizes  noise  effects  due  to  the  high  frequency 
enhancement  which  overcome  MTF  degradation.  The  smoothness,  however,  causes  limited 
performance  on  sharp  edges  restorations. 

We  assume  that  H{u,v)  is  the  average  atmospheric  MTF,  the  original  picture  is  of 

dimensions  MxN ,  and  the  MTF  and  PSF  are  of  dimensions  JxK.  Each  matrix  is  padded  with 
zeros  to  obtain  {M  +  2J)x{N  +  2K)  matrix  size.  A  two  dimensional  second  order  derivation 

is  defined  by  the  Laplacian  matrix 


0  1  0 
1  -4  1 
0  1  0_ 


(12) 


The  mathematical  form  of  the  constrained  filter  is  given  by:^'* 
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M(m,v)  = 


\H{u,v)^\ 


{M  +  2J){N  +  2K)(\H{u,v)^\+Ul{u,v)\ 


(13) 


where  L{u,v)  is  the  Fourier  transform  of  [/],  and  X  is  obtained  by  Lagrange  multipliers 
technique,  or  can  be  evaluated  iteratively.  An  example  of  image  restoration,  obtained  over  a 
6.5  km  horizontal  path  length  through  the  atmosphere,  using  the  constrained  filter  is  presented 
in  Fig.  7.  It  is  clear  that  in  the  restored  image  additional  details  are  resolvable  which  can  not  be 
seen  in  the  original  picture,  such  as  the  two  horizontal  bars  at  the  bottom  of  the  picture. 


5.  Modified  Backus-Gilbert  filter 


The  modified  Backus-Gilbert  filter^^  is  a  spatial  domain  filter.  It  is  assumed  that  the 
atmospheric  point  spread  function  is  a  finite  impulse  response  (FIR),  and  so  is  the  restoring 
filter.  Here,  in  the  examples  to  be  shown  it  is  assumed  that  the  PSFs  dimensions  are  15x15 
pixels.  Two  criteria  in  the  restored  image  are  implemented.  These  are  image  sharpness,  in  order 
to  resolve  details  as  fine  as  possible,  and  minimum  noise  variance.  Image  sharpness  is  obtained 
by  narrowing  the  restored  PSF  as  much  as  possible. 

Defining  the  atmospheric  PSF  as  h{x,y)  and  the  restoring  filter  as  m{x,y),  the 

restored  PSF  is  defined  as 

h{x,y)  =  h{x,y)®m{x,  y).  (14) 

In  order  to  achieve  perfect  restoration,  h{x,y)  should  be  as  close  to  an  impulse  as  possible. 

One  possible  measure  of  a  narrow  restored  PSF  is  a  minimized  second  moment,  i.e. 

R  =  E\^jr^h{x,y)dxdy^  (15) 

where  =x^  +  y^  represents  the  coast  of  spread  of  h{x,y). 

The  requirement  of  minimizing  noise  variance  can  be  stated  mathematically  as: 

V  =  £{n^(x,y)}->min  (16) 

where  n{x,y)  =  n{x,y)®m{x,  y)  (17) 

is  the  noise  of  the  restored  image. 
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Since  minimizing  average  spread  R  and  minimizing  noise  variance  V  can  not  be 
achieved  simultaneously,  a  tradeoff  between  resolution  and  noise  must  be  accomplished,  and  a 
compromise  is  by  minimizing  the  weighted  sum; 

F  =  {l-'k)R  +  XV-^mn,  (18) 

where  A,  should  be  selected  according  to  the  restoration  requirements,  blur  versus  noise.  A 
detailed  mathematical  description  of  solving  (18)  is  presented  in  Ref.  19. 

An  example  of  image  restoration,-  for  the  same  scene  as  in  Fig.  3,  using  the  modified 
Backus-Gilbert  filter,  is  presented  in  Fig.  8.  Another  example  is  in  the  thermal  range,  where 
this  filter  yielded  very  good  restoration  results.  This  is  demonstrated  m  Fig.  9  for  the  same 
scene  as  in  Fig.  5. 


6.  Discussion 

A  comparison  between  the  different  restoration  techniques  indicates  that  the  filter  of 
best  performance  is  the  improved  Wiener  filter.  This  is  explained  by  the  use  of  the  atmospheric 
model  as  a  noisy  filter,  due  to  turbulence  jitter  characteristics.  Thus,  when  dealing  with  high 
resolution  visible  and  near  infrared  imaging  through  long  distances,  atmospheric  distortions 
need  no  longer  be  treated  as  a  deterministic  filter  represented  by  average  MTF,  but  can  now  be 
modeled  as  an  inherent  part  of  the  noise  in  the  received  image.  Determining  this  random  part  is 
essential  to  achieving  improved  restored  images,  and  a  method  of  doing  so  is  presented  here. 
This  method  can  be  used  to  help  overcome  the  jitter  characteristic  of  turbulence,  and  is  capable 
of  real  time  image  restoration  via  parallel  processing  transputers.  As  turbulence  effects 
becomes  weaker,  such  as  in  the  thermal  image  example,  the  advantage  of  this  method 
diminishes,  and  standard  Wiener  filter  performance  is  fairly  good.  The  original  picture's  power 
spectral  density,  for  visible  images,  is  best  estimated  using  a  fractal  model.  This  model  yielded 
very  good  results  using  both  the  fractal  estimation  filter  and  the  improved  Wiener  filter.  This 
model  failed  completely  when  used  to  restore  thermal  images.  It  is  expected  that,  due  to 
energy  transfer  considerations,  thermal  images  will  obey  instead  a  Markovian  power  spectral 
density  model.  Thermal  images  restoration  thus  is  a  subject  for  further  investigation,  using  the 
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Markovian  model.  In  cases  where  a  smooth  image,  rather  than  a  sharp  but  noisy  one  is 
required,  the  constrained  filter  can  be  used  to  satisfy  this  condition.  This  is  also  the  case  when 
using  the  modified  Backus-Gilbert  filter,  which  has  also  the  advantage  of  shorter  computation 
time.  In  any  case  where  turbulence  is  not  too  dominant,  the  last  two  filters  yield  impressive 
results  (as  in  the  thermal  range).  However,  as  atmospheric  MTF  variance  increases,  it  is 
strongly  recommended  to  use  the  improved  wiener  filter,  which  indeed  yielded  the  best  results. 

Conclusions 

Restorations  of  visible  and  thermal  images  distorted  by  the  atmosphere  are  presented 
based  upon  atmospheric  MTF.  Such  restoration  techniques  are  fundamental  and  include  both 
spatial  and  spatial  frequency  filters.  In  the  restored  images  it  is  generally  quite  possible  to  see 
detail  almost  as  small  as  possibly  resolvable  by  instrumentation,  i.e.  pixel  size.  In  the  visible 
image  examples,  pixel  size  at  6.5  km  horizontal  distance  is  -7  cm.  Examples  for  finest 
resolvable  details  are  poles  in  a  fence  and  antenna  bars(~  7  cm).  In  the  thermal  images,  pixel 
size  at  2  km  is  ~62  cm,  and  examples  of  finest  resolvable  details  are  windows  whose  size  is 
about  the  same.  This  means  that,  essentially,  all  atmospheric  blur  is  removed  in  the  restoration 
process,  and  the  distant  object  scene  is  observed  as  if  there  were  no  atmosphere.  This  can  only 
be  possible  by  considering  the  atmospheric  MTF  properly,  including  the  significant  role  played 
by  aerosol  MTF.  Such  restoration  does  not  depend  on  target  shape.  It  is  a  fundamental  image 
correction,  which  can  be  followed  by  other  image  processing  techniques  as  required.  In  the 
visible  range  an  improved  Wiener  filter  designed  specifically  for  imaging  through  the 
atmosphere,  together  with  a  fractal  estimation  model  for  the  object's  power  spectral  density, 
appears  to  give  excellent  results.  In  the  thermal  infrared,  where  turbulence  effect  is  less  severe, 
several  different  filters  each  give  excellent  results.  There  are  cases  where  turbulence  is  not  so 
severe,  and  aerosol  MTF  plays  the  most  important  role  in  determining  the  atmospheric  MTF, 
such  as  in  the  case  of  thermal  imaging,  and  the  image  restoration  task  is  much  easier  since  the 
atmospheric  MTF  then  is  nearly  deterministic.  Here,  however,  the  correction  is  for  both 
aerosol  and  turbulence  derived  blur. 
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Figure  captions 

Fig.  1.  Restoration  of  two  scenes  at  6.5  km  horizontal  distance  (visible  wavelengths)  using 
fractal  PSD  best  fit. 

Fig.  2.  Restoration  of  a  thermal  (8-12  pm)  image  of  a  scene  at  2  km  horizontal  distance 
using  fractal  PSD  best  fit. 

Fig.  3.  Improved  Wiener  filter  restoration,  using  fractal  PSD  estimation  of  scenes  of  Fig.  1. 

Fig.  4.  Standard  Wiener  filter  restorations  of  scenes  of  Fig.  1. 

Fig.  5.  Improved  Wiener  filter  restoration  of  a  thermal  (8-12  pm)  image  at  2  km  horizontal 
distance. 

Fig.  6.  Standard  Wiener  filter  restoration  of  the  thermal  image  of  Fig.  5. 

Fig.  7.  Restoration  of  scene  at  6.5  km  horizontal  distance  (visible  wavelengths)  using 
constrained  filter. 

Fig.  8.  Restoration  using  modified  Backus-Gilbert  filter  of  scene  of  Fig.  3. 

Fig.  9.  Restoration  of  a  thermal  image  using  modified  Backus  Gilbert  FIR  filter  (same  scene 
as  Figs.  5  and  6). 
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original  picture  restoration  with  FIR  filter 
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